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Analysis of the exoplanet containing system Kepler-91 
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Abstract We have applied the graphical user in¬ 
terfaced close binary system analysis program Win- 
Fitter to an intensive study of Kepler-91 using all 
the available photometry from the NASA Exoplanet 
Archive (NEA) at the Caltech website: 

http: //exoplanetarchive.ipac.caltech.edu 

Our fitting function for the tidal distortion derives from 
the relevant Radau equation and includes terms up to 
the fifth power of the fractional radius. This results in 
a systematic improvement in the mass ratio estimation 
over that of Lillo-Box et al. (2014a) and our derived 
value for the mass ratio is in close agreement with that 
inferred from recent high-resolution spectroscopic data. 

It is clear that the data analysis in terms of simply an 
eclipsing binary system is compromised by the presence 
of significant other causes of light variation, in partic¬ 
ular non-radial pulsations. We apply a low-frequency 
filtering procedure to separate out some of this addi¬ 
tional light variation. Whilst the derived eccentricity 
appears then reduced, an eccentric effect remains in 
the light curve. We consider how this may be main¬ 
tained in spite of likely frictional effects operating over 
a long time. There are also indications that could be 
associated with Trojan or other period-resonant mass 
concentrations. Suggestions of a possible secular period 
variation are briefly discussed. 


E. Budding 
q;. Puskiillii 
M. D. Rhodes 
O. Demircan 
A. Erdem 

^University of Canakkale, TR 17020, Turkey; 
^Carter Observatory and 

^SCPS, Victoria University of Wellington and 
^Dept. Physics & Astronomy, UoC, New Zealand. 
®BYU, Provo, Utah 


Keywords stars - close binary; stellar dynamics; exo¬ 
planets; light curve analysis 

1 Introduction 

Astrophysical aims of the Kepler Mission, together with 
practical details on its construction were set out by 
Borucki et al. (2003). Devore et al. (2009) provided 
general motivational background, noting the 400 year 
anniversary of Kepler’s publication of the laws of ellip¬ 
tical orbits and equal areas in the year of the satellite’s 
launch. The Ames Research Center has had a promi¬ 
nent role in the materialization of such purposes. A 
comprehensive early summary was that of Borucki et 
al. (2011). In mid-2013 NASA announced that two out 
of the original four reaction wheels used in pointing the 
telescope on the satellite had become inoperable and 
the initial objectives of the mission would be compro¬ 
mised as a result. However, there is, by now, a large 
archive of photometric data inviting continued close at¬ 
tention and discussion. 

Rhodes & Budding (2014), giving further back¬ 
ground information relevant also to the present paper, 
tested their light-curve fitting software for Kepler ex¬ 
oplanet transit light curves against results published 
by others. They concluded that there was a fair mea¬ 
sure of agreement about published parameter values al¬ 
though significant differences were found in some cases. 
Many of these differences can be associated with the 
transition of the photometric to inferred absolute pa¬ 
rameters. Rhodes & Budding (2014) accounted for this 
in terms of the high sensitivity of derived mean densi¬ 
ties to imprecisely known, but observationally obtained, 
surface gravities. In the case of KOI 13.01, however, 
it is feasible that the differences in derived inclina¬ 
tions are valid and can be associated with real 
short term changes due to precession (Szabo et 
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al., 2013). KOI 377.01 (= Kepler-9) is similarly open to 
possible variations of empirical parameters associated 
with orbit complexities of this known multi-planet sys¬ 
tem. For KOI 3.01, the difference in the derived stellar 
relative radius (ri) may again arise from orbital com¬ 
plexity (Bakos et al. 2010). The unusually high limb- 
darkening coefficient found for Kepler-1 has remained 
with further study of the data and is an issue calling 
for continued attention and analysis, perhaps related 
to other peculiarities associated with this hot Jupiter 
(Kipping & Spiegel, 2011). 

The fitting program WinKepler, discussed in 
Rhodes & Budding (2014), has since been upgraded and 
is now called WinFitter. It performs optimization by 
a modified Marquardt-Levenberg application of a fit¬ 
ting function to a photometric data-set (light curve). 
The fitting function is based on the Radau model devel¬ 
oped from Kopal’s (1959) approach to the tidal and 
rotational distortions (ellipticity), together with the 
radiative interactions (reflection), of massive and rel¬ 
atively close gravitating bodies. It can be downloaded 

from http://home.comcast.net/~michael.rhodes/, 

The upgrade allows the user the option of regular 
stellar eclipsing binary light curve analysis as well as 
that of exoplanet light curves such as those released by 
the Kepler Mission via the Mikulski Archive for Space 
Telescopes (MAST). 

Optimal models produced by WinFitter corre¬ 
spond to the least value of defined as S(^o,i — 
lc,i)^/Alf (Bevington, 1969), where lo,i and 4,1 are 
the observed and calculated light levels at a particular 
phase. Ali is an error estimate for the measured val¬ 
ues of lo,i- The NASA Exoplanet Archive (NEA) lists 
empirical values of Ali for each datum, which guide 
the values assigned to WinEitter. Given such start¬ 
ing A4s, we have retained flexibility in the values to 
be used in the calculation, however, as it became 
clear that there are sources of photometric variation 
other than the eclipsing binary effects and photon noise. 
The question of what can be regarded as random noise, 
from the point of view of systematic effects on the close 
eclipsing system model, is a matter we discuss further 
in Section 2.2. 

The theoretical light-level 4 corresponds to the given 
fitting function. The central problem in the analysis 
of photometric data, like that discussed in the present 
paper, is to find the best values of the parameters of 
this fitting function, together with the demonstration 
of a formally determinate optimal set of such parame¬ 
ters, whilst ensuring adequacy of the underling model 
to account for the observed effects. Location of a sin¬ 
gle optimum in the hyperspace formed by the unknown 
parameters and is effected by the simultaneous van¬ 
ishing of the y^ derivatives with respect to all those 


sought parameters. Determinacy is checked by the y^ 
Hessian staying positive definite in the vicinity of this 
optimum. Adequacy of the model entails that the value 
of y^/i^, where v is the number of degrees of freedom 
of the data-set, should be acceptably close to unity at 
the optimum. Standard tabulations of the 
ues at given probability levels can be used to check 
this. Rhodes & Budding (2014) included further dis¬ 
cussion of this subject as well as the physical basis of 
the modelling (cf. also Budding & Demircan, 2007 (Ch. 
9); Rhodes, 2015). 

2 Kepler-91 

2.1 Preliminary information 

Kepler-91 (= KOI 2133.01; KIC 8219268 0 is be¬ 
lieved to contain an unusual combination of a post- 
Main-Sequence star ascending the Red Giant Branch 
accompanied by a jupiter-sized planet heated to around 
2000 K (Huber et al., 2013; Lillo-Box et ah, 2014a). 
Planetary transits are not immediately apparent in in¬ 
dividual long-cadence data-sets: the identification was 
achieved by Tenenbaum et al. (2013) after applying 
special periodicity-finding procedures. We used the 
ephemeris attributed to Tenenbaum et al. (2013), as 
listed in the NEA, in phasing the data and the adopted 
epoch and period are given in Table 1. In doing this 
we have assumed that the adopted Epoch does 
indeed correspond to an identified transit mid¬ 
minimum. Some tens of giant planet -I- giant star 
configurations have been found hitherto and they have 
attracted interest related to the theory of their origin 
and evolution (cf. e.g. Lin et ah, 1996; Johnson et ah, 
2007; Nagasawa et ah, 2008; Villaver & Livio, 2009). 

An aspect of this subject concerns possible engulf- 
ment of a close planet as the host star expands at the 
end of its Main Sequence stage. That Kepler- 9110 is a 
hot jupiter in such a situation was confirmed by Lillo- 
Box et al. (2014a,b), and more recently by Barclay et 


^KIC stands for Kepler Input Catalogue, KIC for Kepler 
Object of Interest. The webpage http://exoplanetarchive 
ipac.caltech.edu/docs/data.html7redirected guides users 
to numerous data-files, including the ‘cumulative list’ 
that we used to obtain prior parameter estimates. 
The official website for Kepler light curves is actually 
https://archive.stsci.edu/kepler/downloads_options.html. How¬ 
ever, the NASA Exoplanet site referred to above does not require 
file-type conversion, and has convenient normalization and useful 
plotting options. The same light curve data can be downloaded 
from either site. 

^The symbol b is associated with the hot jupiter component. 
Kepler-91 may refer to the entire system or sometimes just the 
host star. 
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al. (2015), though that picture had been contested by 
Esteves et al. (2013), and subsequently by Sliski & Kip¬ 
ping (2014). More recently, Esteves et al. (2015) have 
analysed a larger data-set and retracted from their pre¬ 
vious findings. Their latest results are included in Ta¬ 
ble 4 for comparison. Some of the earlier discussion 
related to density estimates, whose sensitivity to er¬ 
ror may be associated with the cubed parameter ri in 
the denominator of the constraint p* = 3/47rGP^ri^. 
If an absolute radius i?* is estimated from separate 
evidence, the density is alternatively constrained by 
p* = 3p/47rGi?*, but again a relatively large propor¬ 
tional error in the observationally determined value of 
the surface gravity g would compromise the density esti¬ 
mate (Muirhead et ah, 2012; Rhodes & Budding, 2014). 

Lillo-Box et al. (2014a), using the 2.2 m telescope 
at Calar Alto (Spain), checked for the possibility of 
background contamination associated with the rela¬ 
tively large angular size of the Kepler pixels (^4 arcsec). 
They also looked for possible non-planet-related trends 
in individual data segments. Commonly used flux mea¬ 
sures have been tabulated by the archive data-source 
under the heading PDCSAP (pre-search data condi¬ 
tioning simple aperture photometry) fluxes. These data 
have resulted from additional processing after the SAP 
data (also tabulated), with the aim of mitigating data 
artefacts. Sometimes, investigators apply separate cor¬ 
rection procedures on the SAP data, but Lillo-Box et 
al. (2014a) deemed it sufficient to utilize the listed PD¬ 
CSAP information. Barclay et al. (2015) also started 
with the PDCSAP fluxes, though they subsequently 
modelled and separated out systematic effects in the 
data. 

In preparing their analysis of the system, Lillo-Box 
et al. (2014a) noted the importance of reliable host star 
parameters, particularly where these can be elucidated 
by separate means prior to study of the Kepler photom¬ 
etry. They made use of high resolution spectroscopy of 
Kepler-91 from Calar Alto to this end, fitting, for exam¬ 
ple, model spectral energy distributions, matching indi¬ 
vidual spectral lines and including detailed asteroseis- 
mological analysis. In this way, they were able to give 
estimates of the radius and mass values as 6.3±0.16 Rq 
and 1.31±0.10 Mq. Other relevant physical attributes 
were given too, including an age, estimated to be not 
far from 5 Gy, though with a substantial margin of er¬ 
ror. Table 2 of Lillo-Box et al. (2014a) collects together 
stellar parameters from a variety of previous sources, 
where a fairly wide range of values can be noted (ap¬ 
preciably larger than the quoted errors). Our Table 1 
assembles input information needed to build up a more 
comprehensive picture of the system. 

The fact of pulsational activity in this star, relating 
to the aforementioned asteroseismological work, could 


be anticipated on the basis of the arguments used 
in Rhodes & Budding (2014) concerning the level of 
scatter observed around the mean value of the stel¬ 
lar flux. Rhodes & Budding (2014) found that for a 
steady source affected only by Poissonian noise, this 
should work out at ^ ^/^/23.8, where the Poisson factor 
ip = Nf being the mean PDCSAP flux count, 

in this case for the long cadence data-sets. This would 
lead to a datum probable error of about 0.00012 (120 
ppm), whereas what is found in typical light curve fit¬ 
ting is a scatter (regarded as noise) of about 3 times this 
level. There are thus clear indications in these data-sets 
of effects that regular eclipsing binary modelling does 
not take into account. The analysis of Barclay et al. 
(2015) addressed this point in some detail, though their 
light-curve rectification to an eclipsing binary model 
was on an essentially empirical basis. 

2.2 Kepler photometry analysis 

2.2.1 Preparation 

In our analysis of this system using WinEitter, 
we used all available long cadence data sets for the 
observing runs (quarters) 1-17, which we downloaded 
from the NEA. For each data-set we processed the 
listed data points using the period of Tenenbaum et 
al. (2013) and normalized the PDCSAP fluxes to unity. 
The given times of observation (BKJD) were then con¬ 
verted to phases in the range 0.0 to 1.0. In this way, 
we constructed working light curves of ^4000 points. 
These data-sets proved informative, although they ap¬ 
pear very scattered to the eye. 

For initial values of the input parameters, we could 
use information from the KOI (Kepler Object of Inter¬ 
est: cumulative list) given in the NEA (see also Batalha 
et ah, 2013). We also take into account the parameters 
listed by Lillo-Box et al. (2014a), whose emphasis on 
prior host-star data was mentioned above (see also Sea- 
ger & Mallen-Ornelas, 2003). Our Table 1 numbers are 
thus rounded averages from these sources. A more de¬ 
tailed discussion about prior possibilities for such input 
information was given by Barclay et al. (2015), whose 
analysis procedure allows for variation of hyperparam¬ 
eters. However, the main parameters connecting with 
photometric effects, such as the eclipse and its shape, 
constitute separate information to absolute parameter 
specification. Fittings to the transit minima led to rea¬ 
sonably consistent indications about ri, k and i (the ra¬ 
dius of the star expressed as a fraction of the semi-major 
axis of the orbit, the ratio of planet to stellar radii, and 
the orbital inclination) independently of much of the in¬ 
formation set out in Table 1. The inclination turns out 
to be relatively low (~70° in Table 2) and the star is a 
relatively large fraction (~0.4, Table 2) of the planet- 
star separation, so that a noticeable boundary correc¬ 
tion (cf. Kopal, 1959; ch. 4) can occur for moderate 
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levels of distortion. This can delay the onset of the 
transit by a degree or two of phase from what would 
correspond to a circular stellar disk. The oft-cited for¬ 
mulae of Mandel and Agol (2002), used in transit anal¬ 
yses such as those mentioned above, could thus render 
the fitting function less than adequate for precise data, 
regardless of the precision of available priors. 

Our general procedure is to approach the optimal 
photometric parameter specification in steps. First we 
concentrate on the parameters that have a relatively 
strong effect or are highly constrained by the light 
curve’s form. When these settle towards well-defined 
numbers, we allow for simultaneous improvements in 
weaker parameters. In preliminary fittings we would 
generally start by optimizing the scaling constant for 
the light curve’s vertical axis (unit of light) C/H and 
similarly the x-axis zero point (epoch of eclipse mid¬ 
minimum), or correction to the assigned phases AfjjQ. 
We then use such values in subsequent runs, keeping 
the generally independent A(po fixed, and optimizing 
for [/, ci, k; and i. 

In Figure 1 we show our optimal fitting to the transit 
region. We extracted the transit portions of the light 
curves corresponding to the phase range -0.05 to 0.05. 
Figure I shows the result of combining all the transit 
regions of all 17 light curves into a single sorted data¬ 
set of 6322 points, which has been then binned to a 
representative 97 points in bins of phase interval 0.001. 
Lillo-Box et al. (2014a) report a procedure of initial 
transit analysis: their Fig 6 may be compared with our 
Fig 1. 

2.2.2 More detailed examination 

Using the averaged values for the main geometric 
quantities from the initial fittings (bottom row of Table 
2), we proceeded to examine the complete light curves. 
It can be readily expected from the preliminary infor¬ 
mation that proximity effects due to the relatively low 
separation of this hot Jupiter should be detectable, al¬ 
though the system is fairly faint {V = 12.884, Everett 
et ah, 2012) and the scatter of the data points rela¬ 
tively high, as noted above. In subsequent fittings to 
complete light curves we experimented by concentrat¬ 
ing all the determinability of the data to optimizing 


®In principle, the adopted value of U could be used to scale the 
representative flux from the system and thence derive a corre¬ 
sponding stellar magnitude. However, probably related to the 
wide passband of the Kepler filter, the calibration in question 
turns out not simple, and we can see an appreciable departure 
from the trend of magnitudes with effective wavelength listed in 
Table 3 of Lillo-Box et al (2014a) in the case of the Kepler magni¬ 
tude. Our flux-derived Kepler magnitude (12.408) thus appears 
too faint at the flux-weighted wavelength 0.735 /rm, as judged by 
the run of values cited by Lillo-Box et al. (2014a). 


only U and the mass-ratio, q, the other parameters hav¬ 
ing been fixed from previous fittings. This procedure is 
sometimes called a photometric q-search. The essential 
form of the fitting function used here is as given in Eqn 
(9.17) in Budding & Demircan (2007), with further de¬ 
tails being given in the bibliographical notes (9.7) of 
that book. WinFitter has now multiplied the rel¬ 
ative luminosities Li _2 of the original formulation by 
‘Doppler-beaming’ factors, as explained in the paper of 
Shporer et al (2012). 

We continued with further fittings for the 17 light- 
curves folded over each quarter, allowing the geomet¬ 
ric parameters to relax from previously found values. 
Table 2 summarizes the findings of a large number of 
fitting experiments as performed by our team-members 
separately. The run of values through the tables allows 
an insight into the determinability of the parameters, 
apart from the formal errors of each fitting. The 17th 
data set contains only about 1/4 of typical quarters 
and its best-fitting values diverged somewhat from the 
trend of the other solutions, while the reduced re¬ 
mained relatively high. We have therefore included the 
results from that fitting in the tabulated mean values 
with a correspondingly reduced weight. At this point, 
we had failed to confirm any eccentricity to the orbit 
from these separate raw data-sets through finding a sig¬ 
nificant improvement in from inclusion of optimized 
e and ui values. 

We also looked for any trend that might be dis¬ 
cernible in the values of A/jq, which relates to the time 
of transit variation (TTV) exhibited in some exoplanet 
systems. Our results on this are indicated in Figure 2 
with the corresponding numbers in column 3 of Table 
The inclusion of a linear correlation (incorrectly 
assigned period) improves the value of but only 
marginally, and, in fact, x^jv is not significantly al¬ 
tered. A more noticeable reduction in x^ jv (19%) is 
obtained by the parabolic fit (decreasing period) but 
the decrease is from around a one in three ratio of ex¬ 
ceeding the higher ^ jv by pure chance, to around a 
one in three ratio of being less than the lower /v 
value randomly, and so cannot be considered signifi¬ 
cant. Alternatively, the Bayesian information criterion 
(BIG), which is easily related to the variate using the 
formula of Kass & Raftery (1995), gives a very small 
preference for simply keeping the adopted period with 
a slight change of the reference epoch. But the hrst 
3 fittings (reference point only, linear, and parabolic 
trends, have no significant difference in corresponding 


^It should be noted that the effects shown in Fig. 2 are not tim¬ 
ing variations of individual transits, but an average result for 
each quarter. Any trend that might be discerned would not be 
commensurate with the orbital period therefore. 
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BIC values (23.3, 25.0, 24.2; respectively), though a 4- 
parameter fit (BIC = 27.0) would not be supported, by 
the assessment table of Kass & Raftery. 

In fact, the implications of such a period reduction 
would be the highly unexpected circumstance of catch¬ 
ing a planet within a few thousand of years before 
its demise. On the contrary, Sato et al (2015) gave 
a slightly increased period value from that given by 
Tenenbaum et al. (2013) (see Table 4). There are, how¬ 
ever, interesting alternative physical scenarios relating 
to period changes that we discuss later. 

We next combined all ~65000 data points into a rep¬ 
resentative binned light curve that we fitted in various 
ways, starting by using the mean parameters given at 
the bottom of Table 2. The results are summarized in 
Table 3. We can note here the good measure of agree¬ 
ment on the main geometric parameters (ri, k and *), 
as well as the mass ratio, from the various fitting ap¬ 
proaches. 

Lillo-Box et al. (2014a) carried out a basically simi¬ 
lar procedure, aiming also to determine a planet to star 
mass ratio, and thence, using their separate information 
on the star, a planetary mass. Their fitting function 
combines the three contributions of ellipticity, Doppler 
beaming (a relatively very small term) and reflection, 
given by their Eqns 7, 8 & 9, respectively. Eqns 7 and 
9 are very simple, however,and although their key el¬ 
lipticity contribution derives originally from the same 
Radau-model formulae that we use (Kopal, 1959), the 
neglect of the third and fourth order terms in the tidal 
distortion (given the large relative size of the star) turns 
out to lead to a significant overestimate of the mass 
ratio as determined from the photometric effect of the 
stellar tide. Given that ri « 0.4 we can see that the ne¬ 
glect of the third and fourth order terms will bring their 
mass ratio (0.00064) down to a value in close agreement 
with ours. 

Lillo-Box et al. (2014a) gave the planet’s mass to 
be about 0.88 -I-0.17/-0.33 Mjup but with radius about 
1.384 -l-0.il/-0.054Rjup, i.e. a density of about 1/3 that 
of Jupiter. This was interpreted as an atmospheric in¬ 
flation, associated with the relatively high stellar irra¬ 
diation. They also reported an eccentric orbit (e = 
0.066 -I-0.013/-0.017, w = 316.8-I-21/-7.4 deg), asso¬ 
ciated primarily with an apparent asymmetry of the 
light curve as a whole. They could not confirm this 
definitely from asymmetry of the transit alone, accept¬ 
ing that the eclipse effect was inconclusive regarding 
eccentricity. Our findings are similar on this point (Ta¬ 
ble 3). The complete light curve they modelled consists 
of about 260 individual points that have been binned 
from the original complete data-set, from which they 


had clipped outliers so as to arrive at a working sam¬ 
ple of ~52000 points. From our fittings to the original 
data-sets we estimated the s.d. scatter of such data at 
~0.00037, hence the error-bars shown on the light curve 
of Lillo-Box et al. (2014a) should be typically around 
0.00002 on a Poissonian assumption of their distribu¬ 
tion. This corresponds to the error bars shown in their 
Fig 7. 

It is clear from that Fig 7, however, that the residuals 
are not distributed randomly around the model curve 
with this precision. The model curve goes through quite 
fewer than ~2/3 of the shown error bars, that one could 
reasonably expect it would do on the basis that the data 
corresponds only to the close binary model plus random 
noise. But this would happen anyway as a result of 
the pulsational effects. This matter was addressed by 
Barclay et al. (2015), who introduced an additional em¬ 
pirical noise model that they attributed to the effects 
of surface granulation of the host star. We digitized 
the data points of Fig 7 in Lillo-Box et al. (2014a) and 
could mostly confirm their model parameters to within 
reasonable errors including a similar eccentricity result. 
In order to probe the situation further we binned the 
source data and applied WinFitter to the resulting 
light curve. 

A periodicity analysis was then carried out on the 
residuals from this fitting. We show the result in Figure 
3 (see also the residuals in Figs 4 and 5). A very signifi¬ 
cant peak at ~6/P frequency can be seen, with a lower 
but wider peak at the 3/P submultiple. The corre¬ 
sponding photometric contribution, persisting through¬ 
out the ^240 orbits of the complete data-set, points 
to an orbit-oscillation resonance. This may well cor¬ 
respond to the beat which is produced at ^11.0 /iHz 
between the higher frequency L2 and LO modes shown 
by Lillo-Box et al (2014a), particularly the strong ones 
(n9L02 and nllLOO) between ~ 104 and 115 ^Hz. The 
existence of this additional photometric effect is, of 
course, of interest in its own right; its continuation evi¬ 
dent in the orbitally phased and binned light curve even 
to the eye. The nature of this resonance invites further 
physical discussion that we raise in Section 3. From 
the more immediate aim of parametrization of the close 
binary (star -|- planet) configuration, however, this is a 
complication that it would be desirable to separate out. 
Its retention in the analysis of Lillo-Box et al (2014a) 
may compromise the deduction of an eccentricity effect, 
or the parametrization in general, to some extent. 

We therefore carried out a cleaning operation by fit¬ 
ting the residuals from the original light curve fit with 
a Fourier decomposition reaching to terms in frequency 
6/P. This fitting is shown in Figure 4. Our procedure is 
similar to that used in the cleaning of close binary sys¬ 
tems showing maculation effects (Zeilik et ah, 1988). In 
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Table 1 Primary Input Data 


KOI 

V 

M* 

R* 

T. 

Epoch 

P 

Z 

logs 

a 

Mp/Mt. 

U 

2133.01 

12.884 

1.3 

6.4 

4600 

2454969.3966 

6.246580 

0.11 

2.94 

0.073 

0.00005 

0.74 


We have here used the notation: V — conventional V magnitude, M* - mass of star (solar masses), P* - radius of star (solar radii), T* 
- temperature of star (K), P - orbital period (in days), Z - metallicity of star, log g - logio of the surface gravity of star, a — semi-major 
axis in AUs, M^p/M^ — ratio of planet to star masses, u — stellar (linear) limb-darkening coefficient. The numbers given in this table 
are rounded averages from the results of Lillo-Box (2014a) and various sources provided by the NEA, including the cumulative list, 
the confirmed planets and exoplanet transit survey data. 


Table 2 Results of optimal fittings to light curves drawn from the 17 available quarters. 


Qtr 

U 

A 00 

ri 

k 

r2 

i(deg) 

Q 

1 

1.00079 

-1.5 

0.404 

0.0261 

0.0105 

68.2 

0.00039 

2 

1.00084 

-0.5 

0.406 

0.0247 

0.0100 

68.5 

0.00040 

3 

0.99972 

0.2 

0.403 

0.0209 

0.0084 

69.1 

0.00037 

4 

0.99987 

0.3 

0.403 

0.0206 

0.0083 

70.0 

0.00021 

5 

0.99985 

-0.6 

0.403 

0.0229 

0.0092 

70.1 

0.00039 

6 

0.99991 

2.0 

0.407 

0.0228 

0.0092 

68.6 

0.00065 

7 

0.99980 

0.6 

0.403 

0.0183 

0.0074 

69.8 

0.00048 

8 

0.99981 

0.6 

0.400 

0.0224 

0.0089 

70.4 

0.00046 

9 

0.99984 

-0.4 

0.393 

0.0215 

0.0084 

70.3 

0.00027 

10 

0.99984 

0.7 

0.406 

0.0229 

0.0093 

69.1 

0.00051 

11 

0.99988 

-2.1 

0.405 

0.0237 

0.0096 

68.9 

0.00054 

12 

0.99984 

0.5 

0.396 

0.0220 

0.0087 

69.5 

0.00042 

13 

1.00000 

-0.6 

0.395 

0.0220 

0.0087 

70.1 

0.00012 

14 

1.00002 

-0.9 

0.403 

0.0221 

0.0089 

69.6 

0.00036 

15 

0.99994 

0.6 

0.405 

0.0206 

0.0084 

68.8 

0.00056 

16 

0.99967 

-1.7 

0.407 

0.0213 

0.0086 

69.9 

0.00048 

17 

0.99948 

-0.5 

0.409 

0.0196 

0.0080 

69.2 

0.00055 

Mean 

0.99999(34) 

-0.2(1.0) 

0.4026(43) 

0.0220(18) 

0.0089(7) 

69.4(7) 

0.00045(6) 
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effect, we treat the light-curve oscillatory disturbance 
as a low-frequency superposition. The higher frequency 
components that may also be present become strong at 
'^10 times our frequency limit, but such high frequency 
components are more akin to noise in their effects on 
the eclipsing binary model. Our final cleaned model 
and its fitting are shown in Figure 5. We may note that 
the cleaned light curve shows a reduced eccentricity: a 
direct fitting of the raw data may thus be be open to 
doubt if interpreted simply as a real elliptical orbit. At 
least, its determination in a data-set containing both 
eclipsing binary and other causes of photometric vari¬ 
ability cannot be taken at face value. 

In Table 4 our final results and their implications for 
absolute parameters are listed and compared with the 
results of other authors. 


3 Discussion 

3.1 Mass ratio 

The preceding presentation on the Kepler Mission’s 
photometric data on the Kepler-91 system to a large 
extent supports the picture given by Lillo-Box et al. 
(2014a, b), Barclay et al. (2015), Sato et al. (2015) and 
Esteves et al. (2015). The close binary system pho¬ 
tometry analysis program WinFitter that we have 
adapted and applied, however, includes a fuller descrip¬ 
tion of the surface distortions due to tides and rotation, 
as well as the reflection effect, than has generally been 
used in similar studies hitherto. Given the relatively 
large fractional size of the star in the present case the 
higher power terms are significant and we can expect a 
systematic improvement in the mass ratio (reduced by 
^ 40% from that determined by Lillo-Box et ah, 2014a). 
Our derived value q ^ 0.00044 ± 0.00006 for the mass 
ratio is then in good agreement with that which can be 
inferred from the relatively high-accuracy radial veloc¬ 
ity data of Barclay et al. (2015) and perhaps even more 
with the greater coverage and higher mean accuracy 
results of Sato et al. (2015). 

3.2 Additional effects 

It is clear from the comparison of the observed scat¬ 
ter in the data compared with the photon noise ex¬ 
pectable for a K « 12.88 star that there are additional 
short term variations, which can be associated with the 
non-radial pulsational spectrum also studied by Lillo- 
Box et al (2014a). In fact, the scatter in the raw data 
makes individual planetary transits only Itr events: bin¬ 
ning by a factor of ~100 therefore appears desirable to 
allow a reasonably confident approach to parametriza- 
tion. From the empirical point of view, the light curve 
produced in this way can be characterized by, in addi¬ 
tion to regular eclipsing binary system effects, (a) an 
‘O’Connell effect’ (asymmetric maxima) and (b) the 
persistence of quasi-periodic small light drops or dims. 
The former can be associated with the effect of orbital 
eccentricity, given the slight, though calculable, effects 
of a correspondingly varying tidal distortion of the star. 
The possibility of photometric effects associated with 
some surface maculation may be mentioned, given the 
cool temperature and likely structure of this star. We 
know, however, from the slow rotation speed derived by 
Lillo-Box et al (2014a) that the planet performs around 
7.6 revolutions in one rotation of the star (unlike the 
synchronized arrangement usually found in close binary 
stars - cf. e.g. Zahn, 1977). Maculation seems thus un¬ 
likely to persist as coherent systematic effects in the 
binned light curve in general, though starspots may be 
partly associated with the relatively high scatter. 



X^/i^ A/ E5 


Table 3 Results of optimal fittings to various light curve data-sets. 


A0O El k E2 r 2 E3 i(deg) q E4 e cu 


TO 

-0.15(4) 

4.08(7) 

2.21(3) 

8.9(2) 

68.9(5) 






BF 

-0.2(3) 

4.03(6) 

2.16(2) 

8.7(2) 

69.6(5) 

4.2(8) 



1.19 

2.5 

BFe 

-0.3(2) 

4.09(4) 

2.22(2) 

9.04(8) 

69.8(4) 

4.4(6) 

0.05(2) 

300(15) 

1.07 

2.5 

BTO 

-0.1(1) 

4.1 (2) 

2.2(1) 

8.9(4) 

68.7(8) 

7(7) 



1.01 

2.2 

BTOe 

-0.1 

4.1 

2.1 

8.8 

69 

6 

0.03 

260 

1.01 

2.2 

NEAK 


3.8(2) 

2.159(+8,-31) 

8.2 

71.05 






NEAC 





67.4(6) 

4.8(6) 





LB 


4.08(+l,-2) 

2.26(+3.-10) 


69(1) 

6.4(+l,-2) 






Table 3 notation 

TO = Transit Only: means from the individual runs of each quarter 

BE = Binned Full: fitting to a single binning down to 360 points 

BFe = Binned Full: as above, but with eccentricity allowed 

BTO = Binned Transit Only: fitting of all transits, binned to 256 points 

BTOe = Binned Transit Only: as above, but with eccentricity allowed 

NEAK = NASA Exoplanet Archive KOI Cumulative list 

NEAC = NASA Exoplanet Archive - Confirmed Planets 

LB = Lillo-Box et ah, (2013) 

Formal error estimates, where available, are indicated by the digits in parentheses. For reasons of space and convenient comparison 
in one panel, these are here truncated to the rightmost digits in the parameter values. En means that the listed value is 10^ times 
the actual value. The unit of light values {U) for these fittings (not given here) are close to unity, with the very low scatter consistent 
with the errors indicated in Table 3 (i.e. ~23 ppm for Al or ~1 ppm for AC/). 

The BFe row of this table is adopted for our final presentation in Table 4. We present the reduced x^/^ value for this and the BE 
solution with the same error estimate of 0.000025 for the normalized 360 point binned data-set. The decrease of jv from 1.19 to 
1.07 is significant at the 95% confidence level, i.e. there is a relatively high probability to to the improvement produced by an eccentric 
orbit model. This is not the case for the transit-only fittings, where, though modelling is effective in reducing the relative scatter of 
the residuals, the parameters are less well defined. In fact, if eccentricity is allowed to be simultaneously adjustable the x^ Hessian 
becomes non-positive-definite, so that regular error estimates lose meaning. 


Table 4 We list here the Lillo-Box et al. (2014a) (LB) parameters, as well as those of Barclay et al. (2015), Sato et al. 
(2015) and Esteves et al. (2015), alongside those of the present paper (P). Error estimates, influenced by the compilation of 
results, are listed in the final column. Where no error is listed the corresponding parameter has been taken from information 
separate to the fitting. 


Parameter 

LB 

Barclay 

Sato 

Esteves 

P 

Error 

Epoch 

2454969.3966 

2454969.3837 

2454969.3866 

2454969.3958 

— 


P (d) 

6.24658 


6.24668005 

6.24658 

— 


dist. (pc) 

1030 




— 


M, (©) 

1.31 

1.31 

1.31 

1.31 

1.3 


(K) 

2000 


1920 - 2460 

2040 - 3000 

2000 


a 0 

15.62 

15.53 

15.93 

15.71 

15.57 


log gs 

2.95 



2.953 

2.96 


Mp/Ms 

0.00061 

0.00053 

0.00048 

0.00059 

0.00044 

0.00006 

Ts (mean) 

0.403 

0.406 

0.406 

0.401 

0.409 

0.004 

Tp (mean) 

0.0091 

0.0088 

0.0090 

0.0087 

0.0090 

0.0002 

i (deg) 

67-78 

69.17 

67.37 

69.7 

69.8 

0.5 

e 

0.07 

0.02 

0.04 

0.0 

0.05 

0.02 

w (deg) 

320.0 

41.9 

302.8 


300 

20 

b 

0.88 

0.87 

0.90 

0.87 

0.85 

0.03 


In the P column — indicates values as in LB. 
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3.3 Explaining the dims 

Regarding the dims, Lillo-Box et al. raised possible 
lines of explanation involving other planets or satellites 
that co-rotate synchronously or are in a resonant ar¬ 
rangement with Kepler-91b. The idea of eclipses of the 
planet’s reflected light, although apparently modelled 
by Barclay et al. (2015), faces at least conceptual diffi¬ 
culties, since the dims are more than twice as deep as 
the regular scale of reflected light at full phase. When 
this is eclipsed out at the secondary minimum the light 
level should not be appreciably lower than at phases 
just outside the primary minimum; the star’s light only 
being received in either case. We have to countenance 
also improbability going with the paired circumstances 
of another planet’s non-coplanar orbit, though still near 
the line of sight, which is also at a resonant period. 

Given the v-shaped effect at phase ~60°the possibil¬ 
ity of a large Trojan concentration can be considered: 
similarly, perhaps a Hilda concentration could be as¬ 
sociated with the extra light loss about the secondary 
minimum. Such situations have been examined, for ex¬ 
ample by Schwarz et al. (2007), or in the volume edited 
by Souchay & Dvorak (2010), but more detailed clari¬ 
fying explanation is needed. It seems clear that addi¬ 
tional systematic effects are present in the light curve 
that survive the repeated foldings and should there¬ 
fore have some connection to the orbit. However, it is 
also seen that the depth of the dims is ^1/10 that of 
the transit, entailing that any planetoid concentration 
must involve the equivalent of at least several thousand 
Vesta-sized objects if it is to fully account for the scale 
of light loss. The effect is several orders of magnitude 
greater than could apply to the Trojan concentration 
in the Solar System (Yoshida & Nakamura, 2005), and 
so must be regarded cautiously. 

Lillo-Box et al. (2014a) mentioned subtle effects as¬ 
sociated with the Kepler data-processing. A known ef¬ 
fect associated with the relatively large pixel size of 
the Kepler detectors is the presence of a remote Algol 
system in the field. There would again be a paired co¬ 
incidence required in such a case associated with the 
orbital period integral sub-multiplicity. Additionally, 
false-positive scenarios of this kind diminish to low 
probability following high-resolution analysis of data 
from the AstraLux lucky-imaging camera on the 2.2 
m telescope at the Calar Alto Observatory included in 
the study of Lillo-Box et al. (2014a). The more direct 
line of interpretation raised in the previous section is 
thus favoured, though the suggested pipeline error pos¬ 
sibility cannot be entirely ruled out. 

Our treatment, introducing a low pass filter to re¬ 
move these apparently systematic and orbitally res¬ 
onant effects, leads to more precise eclipsing binary 


model defined parameters, as with the approach of Bar¬ 
clay et al. (2015), as well as casting doubt on the pre¬ 
cise representation of light curve asymmetry in terms of 
simply an elliptical orbit. Nevertheless, the light curve 
asymmetry is maintained even after cleaning and the 
most direct explanation of that is in terms of a slightly 
eccentric orbit, even given an expectable trend towards 
orbit circularization (Zahn, 1977), keeping in mind the 
great age of the system estimated by Lillo-Box et al. 
(2014a). At this point, we may note the possibility of 
orbital eccentricity ~ resonant pulsational interactions, 
as considered, for example, in heartbeat’ pulsators (as 
an extreme case: see e.g. Hambleton et ah, 2013). Phys¬ 
ically, the argument is that an orbital eccentricity that 
might otherwise be frictionally eroded can be trapped 
at a certain value if appropriate amounts of energy are 
fed back into the orbit at just the right times from a 
suitable K-mechanism (Henrard, 1982; Alexander, 1988; 
Witte and Savonije, 1997). However, the potentially 
resonant process here is actually a beat effect, rather 
than a direct eigenfrequency of the underlying kappa- 
mechanism. Even so, the amplitude envelope of se¬ 
lected strong eigenfunction pairs does apparently res¬ 
onate at 6/P with the orbital frequency. That the sug¬ 
gested mechanism is, or is not, possible physically, thus 
invites fuller follow-up theoretical consideration. 


3.4 Period variation 

Our modelling suggested the possibility of a steady 
decline in the period, although the likelihood of this is 
not high, given the uncertainties of the time of mini¬ 
mum determination. On the other hand, there appears 
some suggestion, from the results of Sato et al. (2015) 
that the period may be increasing in the short term. 
In a similar way, an apparent trend to decrease of the 
longitude of periastron over the full range of the Kepler 
Mission surveillance cannot be taken so seriously given 
the relatively high errors of the determination. 

But, in connection with such effects, Budding (1984) 
gave a possible form for the variation of the separation 
A between the components of a binary system with to¬ 
tal mass M and two stellar masses mi, m 2 , in which 
mass loss and transfer was taking place as 


A 

A 


2 



mi 

m2 


mi 

m2 


M (e'l -I- £2 -I- £3) 


( 1 ) 


where J is the system’s angular momentum, and the 
small components in Ci are associated with the rota¬ 
tional changes of the two bodies and some surrounding 
matter. These latter components are often neglected, 
at least in a preliminary assessment. 

In classical Algols, regarded as involved in a process 
of mass transfer and where mi denotes the originally 
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more massive and mass-losing star, this equation can 
be simplified to 


^ o [ j {M-2mi) \ 

A I J mi' {M — mi) j 

where we have written simply rn 2 = —mi. In a mass 
loss process, mi is, by definition negative, so this equa¬ 
tion tells us that once mi has dropped below half the 
system’s mass, further mass loss will widen out the sep¬ 
aration, unless there is some significant loss of system 
angular momentum through J/J. This is usually con¬ 
sidered small for classical Algols in the process of con¬ 
servative mass transfer, so the separation of the compo¬ 
nents in this phase is expected to increase. If mi were 
to become very small, but still engaging in mass trans¬ 
fer in this way Eqn 2 has a small denominator, which 
would tend to amplify the separative effect of even a 
low rate of mass transfer. It is feasible that the planet 
Kepler-91b might be losing mass in this way through 
the atmospheric inflation process considered by Lillo- 
Box et al. (2014a). The mean density of the planet 
certainly seems low in comparison to more typical hot 
jupiters. This might then be consistent with the appar¬ 
ent period increase indicated by Sato et al (2015). 

Various possibilities that might produce a noticeable 
variation of J/ J could be considered, including tidal 
friction and magnetic braking. The latter offers a wide 
range of possible parameter variation and appears able 
to introduce significant effects in a variety of situations. 
In general, of course, we should expect such frictional 
effects to make for a negative T/T, which would entail 
period decrease, as suggested in Figure 2. Reasonable 
estimates for J/J due to magnetic braking were esti¬ 
mated by Demircan et al. (2006) and also Erdem & 
Oztiirk (2014). However, the scale of effect that might 
be interpreted from Figure 2 has P jP ~ 10“^, which is 
a 2-3 orders of magnitude greater than the effects con¬ 
sidered in these papers. Granted that the scale of error 
in Fig. 2 makes any present attempt at interpretation 
in physical terms unrealistic, still the possibility of pe¬ 
riod changes is intriguing. A magnetically driven mass 
loss, tending to decrease the orbital period, could be 
reasonably expected from this star; while mass transfer 
from the inflated planet would increase the period. The 
combination of effects in the Kepler-91 system overall 
thus make for interesting alternatives, and the case for 
continued surveillance is compelling. 

This research has been supported by TUBITAK (Sci- 
entihc and Technological Research Council of Turkey) 
under Grant No. 113F353. The help of colleagues in 
the Physics Dept., COMU, including Dr M. Tuysiiz, is 
acknowledged. 



Fig. 1 Binned data covering the transit region for Kepler- 
91b, matched by the WinFitter model. This may be com¬ 
pared with Fig. 6 in Lillo-Box et al. (2014a), or Fig. 4 in 
Esteves et al. (2015). 
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Fig. 2 Differences in observed and determined times of 
minimum for Kepler-91. The straight line (dashed) fit is 
suggestive that the assigned period may have been too long 
- by about 0.00001 d - but the resulting improvement is not 
significant. Somewhat greater significance (in a sense) 
attaches to the secular period decrease (parabola), but even 
there the improvement in probability is not high (see text). 
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Lomb-Scargle Periodicity Analysis 



Fig. 3 Lomb-Scargle analysis of the power in the Kepler- 
91b light curve residuals against orbital frequency. Resid¬ 
uals here come from the differences between the eclipsing 
binary model fitting to the whole, or full-phase, data-set; 
originally of ~65000 points, but binned for this analysis to 
a representative 360 points. The high peak at 6/P is very 
distinct. 



0.00010 4 - 
0.00012 . 



Fig. 5 Complete cleaned and binned light curve of Kepler- 
91, showing proximity effects and an apparent secondary 
minimum, together with the WinFitter model. The resid¬ 
uals are shown below. This may be compared with Fig 7 in 
Lillo-Box et al. (2014a). 


Fig. 4 Fit of a Fourier series with terms up to 6 x the 
orbital frequency. The dims at phases ~l/6 and ~l/2 look 
distinct and may suggest Trojan or Hilda type phenomena. 
This model was added back with a negative sign to the 
original data-set. 
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